MODULE CPTAVE_MOD
CONTAINS
SUBROUTINE CPTAVE(PTHESAT,PTHECAP,PTHEPWP,PCKSAT,PSISAT,PB,PRBARPWP,PRBARCAP)
 
USE PARKIND1  ,ONLY : JPIM     ,JPRB
USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK, JPHOOK

USE YOS_DIM   , ONLY : JPSOILTY, JPTEXT

! (C) Copyright 2005- ECMWF.
!
! This software is licensed under the terms of the Apache Licence Version 2.0
! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
! In applying this licence, ECMWF does not waive the privileges and immunities
! granted to it by virtue of its status as an intergovernmental organisation
! nor does it submit to any jurisdiction.
!     Purpose.
!     --------
!**** *CPTAVE * - COMPUTES AVERAGE SOIL PROPERTIES 
!                 FOR THREE BROAD CLASSES OF SOIL
!           COARSE
!           MEDIUM
!           FINE

!     THE PROPERTIES TO BE DERIVED FOR EACH OF THE THREE SOIL TYPES ARE:
!           SOIL WETNESS AT SATURATION                  (PTHESAT)
!           SOIL WETNESS AT FIELD CAPACITY              (PTHECAP)
!           SOIL WETNESS AT PERMANENT WILTING POINT     (PTHEPWP)
!           CONDUCTIVITY AT SATURATION                  (PCKSAT)
!           SOIL WATER POTENTIAL AT SATURATION          (PSISAT)
!           B EXPONENT FOR CH FORMULA                   (PB)

!     ASSUMPTIONS:
!        1) CLAPP AND HORNBERGER (1978, CH) TYPE OF DEPENDENCY OF PSI AND
!           K ON SOIL WETNESS
!        2) SATURATION SOIL WETNESSES (POROSITY) FROM CH
!        3) PWP AND CAP SOIL WETNESSES FROM PATTERSON (1990, PAT)

!     METHOD:
!        1) GET ZWSAT FROM CH, ZWCAP AND ZWPWP FROM PAT FOR THE 11
!           USDA TEXTURAL CLASSES
!        2) AVERAGE (1) VALUES FOR EACH OF THE 3 CLASSES
!              COARSE        1,2,3,6
!              MEDIUM        4,5,7,8
!              FINE          9,10,11
!        3) GET ZKSAT FROM CH FOR EACH OF THE TEXTURAL CLASSES
!        4) AVERAGE (3) FOR THE CLASSES DEFINED IN (2)
!        5) DETERMINE PSISAT, B FOR EACH OF THE 3 CLASSES BY SOLVING
!           SIMULTANEOUSLY
!              PSI(PTHEPWP)=-15 BAR
!              PSI(PTHECAP)=-0.33 BAR
!        Explicit arguments :
!        --------------------

!        Implicit arguments :
!        --------------------
!     Externals.
!     ----------
!         AVECLA
!     Author.
!     -------

!     Modifications.
!     --------------
!        Original :
!        MODIFICATION  : J.F. Estrade *ECMWF* 03-10-01 move in surf vob
!        M.Hamrud      01-Oct-2003 CY28 Cleaning

IMPLICIT NONE

!     CONSTANTS FOR DEFINING THE SOIL WATER HOLDING CHARACTERISTICS

!     THE FOLLOWING ARE FUNCTION OF SOIL TYPE

!        ZWSAT - CLAPP AND HORNBERGER, 1978
!        ZWCAP, ZWPWP - PATTERSON, 1990

!           (NOTE REVERSAL OF POSTION 4 AND 5 IN PATTERSON TABLES)

REAL(KIND=JPRB) , DIMENSION(JPSOILTY) :: ZWSAT=(/&
 & 0.395_JPRB, 0.410_JPRB, 0.435_JPRB, 0.485_JPRB, 0.451_JPRB,&
 & 0.420_JPRB, 0.477_JPRB, 0.476_JPRB, 0.426_JPRB, 0.492_JPRB,&
 & 0.482_JPRB/)  
REAL(KIND=JPRB) , DIMENSION(JPSOILTY) :: ZWCAP=(/&
 & 0.093_JPRB, 0.172_JPRB, 0.229_JPRB, 0.341_JPRB, 0.288_JPRB,&
 & 0.280_JPRB, 0.345_JPRB, 0.317_JPRB, 0.311_JPRB, 0.338_JPRB,&
 & 0.384_JPRB/)  
REAL(KIND=JPRB) , DIMENSION(JPSOILTY) :: ZWPWP=(/&
 & 0.033_JPRB, 0.058_JPRB, 0.100_JPRB, 0.153_JPRB, 0.145_JPRB,&
 & 0.157_JPRB, 0.201_JPRB, 0.187_JPRB, 0.188_JPRB, 0.207_JPRB,&
 & 0.250_JPRB/)  

!     CONSTANTS FOR HYDRAULIC DIFFUSIVITY AND HYDRAULIC CONDUCTIVITY

!     THE FOLLOWING ARE FUNCTION OF SOIL TYPE (CLAPP AND HORNBERGER)

REAL(KIND=JPRB) , DIMENSION(JPSOILTY) :: ZWCONS=(/&
 & 176.E-6_JPRB, 156.E-6_JPRB, 34.7E-6_JPRB, 7.20E-6_JPRB, 6.95E-6_JPRB,&
 & 6.30E-6_JPRB, 1.70E-6_JPRB, 2.45E-6_JPRB, 2.16E-6_JPRB, 1.03E-6_JPRB,&
 & 1.28E-6_JPRB/)  
REAL(KIND=JPRB) :: ZPSIPWP(JPTEXT),ZPSICAP(JPTEXT)

REAL(KIND=JPRB) , INTENT(  OUT) :: PTHESAT(:)
REAL(KIND=JPRB) , INTENT(  OUT) :: PTHECAP(:)
REAL(KIND=JPRB) , INTENT(  OUT) :: PTHEPWP(:)
REAL(KIND=JPRB) , INTENT(  OUT) :: PCKSAT(:)
REAL(KIND=JPRB) , INTENT(  OUT) :: PSISAT(:)
REAL(KIND=JPRB) , INTENT(  OUT) :: PB(:)

REAL(KIND=JPRB) , INTENT(  IN ) :: PRBARPWP
REAL(KIND=JPRB) , INTENT(  IN ) :: PRBARCAP

INTEGER(KIND=JPIM) :: J

REAL(KIND=JPRB) :: ZFAC
REAL(KIND=JPHOOK) :: ZHOOK_HANDLE

!     THE FOLLOWING ARE THE USDA SOIL TYPES
!     SOIL TYPE 1 (SAND)
!     SOIL TYPE 2 (LOAMY SAND)
!     SOIL TYPE 3 (SANDY LOAM)
!     SOIL TYPE 4 (SILT LOAM)
!     SOIL TYPE 5 (LOAM)
!     SOIL TYPE 6 (SANDY CLAY LOAM)
!     SOIL TYPE 7 (SILTY CLAY LOAM)
!     SOIL TYPE 8 (CLAY LOAM)
!     SOIL TYPE 9 (SANDY CLAY)
!     SOIL TYPE 10 (SILTY CLAY)
!     SOIL TYPE 11 (CLAY)
!     -----------------------------------------------------------------
!*       1.  COMPUTE PTHESAT,PTHECAP,PTHEPWP.
!            --------------------------------
IF (LHOOK) CALL DR_HOOK('CPTAVE_MOD:CPTAVE',0,ZHOOK_HANDLE)
CALL AVECLA(ZWSAT,PTHESAT)
CALL AVECLA(ZWCAP,PTHECAP)
CALL AVECLA(ZWPWP,PTHEPWP)

!     -----------------------------------------------------------------
!*      2.   COMPUTE KSAT.
!            -------------
CALL AVECLA(ZWCONS,PCKSAT)

!     -----------------------------------------------------------------
!*      3.   COMPUTE PSISAT, PB.
!            -------------------
ZFAC=1000._JPRB*100._JPRB/(1000._JPRB*9.8_JPRB)
ZPSIPWP(1)=PRBARPWP*ZFAC
ZPSIPWP(2)=PRBARPWP*ZFAC
ZPSIPWP(3)=PRBARPWP*ZFAC
ZPSICAP(1)=PRBARCAP*ZFAC
ZPSICAP(2)=PRBARCAP*ZFAC
ZPSICAP(3)=PRBARCAP*ZFAC
DO J=1,JPTEXT
  PB(J)=LOG(ZPSICAP(J)/ZPSIPWP(J))/LOG(PTHEPWP(J)/PTHECAP(J))
  PSISAT(J)=ZPSICAP(J)/((PTHESAT(J)/PTHECAP(J))**PB(J))
ENDDO

!     -----------------------------------------------------------------

IF (LHOOK) CALL DR_HOOK('CPTAVE_MOD:CPTAVE',1,ZHOOK_HANDLE)
CONTAINS

SUBROUTINE AVECLA(PSOILAV,PTEXTC)
!     AVERAGES ARRAY PIN(11) USDA SOIL CLASSES INTO
!     ARRAY POUT(3) BROAD TEXTURAL CLASSES
!     SEE PATTERSON (1990)

!     DUMMY ARGUMENTS
REAL(KIND=JPRB),INTENT(IN)  :: PSOILAV(:)
REAL(KIND=JPRB),INTENT(OUT) :: PTEXTC(:)
REAL(KIND=JPHOOK) :: ZHOOK_HANDLE

IF (LHOOK) CALL DR_HOOK('CPTAVE_MOD:AVECLA',0,ZHOOK_HANDLE)
PTEXTC(1)=(PSOILAV(1)+PSOILAV(2)+PSOILAV(3)+PSOILAV(6))/4._JPRB
PTEXTC(2)=(PSOILAV(4)+PSOILAV(5)+PSOILAV(7)+PSOILAV(8))/4._JPRB
PTEXTC(3)=(PSOILAV(9)+PSOILAV(10)+PSOILAV(11)     )/3._JPRB
IF (LHOOK) CALL DR_HOOK('CPTAVE_MOD:AVECLA',1,ZHOOK_HANDLE)

END SUBROUTINE AVECLA

END SUBROUTINE CPTAVE
END MODULE CPTAVE_MOD
